#### prepare packages and functions ####
require(estimatr)
require(dichromat)

issue.var.name <- c("a.article9", "b.defense", "c.revisionism", "d.women", 
                    "e.gay", "f.immigrant", "g.growth", "h.tax")

color <- colorschemes$Categorical.12
plot.col <- c(color[12], color[12], color[12], 
              color[10], color[10], color[10], 
              color[2], color[2])

issue.label <- c("Revision of Article 9", "Increasing Defense Power", 
                 "Historical Revisionism", "Women's Empowerment", 
                 "Gay Marriage", "Accepting Foreign Workers", 
                 "Economic Growth", "Progressive Tax")

# function to compute MMs
MM <- function(data, attribute, outcome, id) {
  lm.data <- data.frame(id = data[, id], A = data[, attribute], 
                        Y = data[, outcome])
  l <- nlevels(lm.data$A)
  result.matrix <- matrix(NA, l, 3)
  for (i in 1:l) {
    level.label <- levels(lm.data$A)[i]
    result.matrix[i, ] <- unlist(lm_robust(Y ~ 1, data = lm.data, 
                                           subset = A == level.label, 
                                           clusters = id)[c(1, 6, 7)])
  }
  result.matrix
}

# function for equivalene testing
equivarence.test <- function(first.group,  # data frame of the first group
                             second.group,  # data frame of the second group
                             covariate,  # covariates' names
                             epsilon  # equivalence limit
) {
  x.t <- first.group[, covariate]  # covariate in the first group
  x.c <- second.group[, covariate]  # covariate in the second group
  n.t <- sum(! is.na(x.t))  # omitting missing values
  n.c <- sum(! is.na(x.c))
  mean.t <- mean(x.t, na.rm = TRUE)  # mean in the treatment group
  mean.c <- mean(x.c, na.rm = TRUE)  # variance in the control group
  var.t <- var(x.t, na.rm = TRUE)  # mean in the treatment group
  var.c <- var(x.c, na.rm = TRUE)  # variance in the control group
  pooled.sd <- sqrt(((n.t - 1) * var.t + (n.c - 1) * var.c) / (n.t + n.c - 2))  # pooled standard deviation
  t.stat <- sqrt(n.t * n.c * (n.t + n.c - 2) / (n.t + n.c)) * (mean.t - mean.c) / 
    sqrt(sum((x.t - mean.t) ^ 2, na.rm = TRUE) + sum((x.c - mean.c) ^ 2, na.rm = TRUE))  # test statistic
  ncp <- (n.t * n.c * epsilon ^ 2) / (n.t + n.c)  # noncentrality parameter
  std.diff <- (mean.t - mean.c) / pooled.sd  # standardized difference in means
  p.values <- pf(abs(t.stat) ^ 2, 1, n.t + n.c - 2, ncp)  # p-value in equivalence testing
  c(std.diff, p.values)
}

# wrapper function for equivalence tests
EQT.wrapper.label <- function(data.tested,  # data frame containing variables
                              treatment,  # treatment variable's name
                              covariates,  # covariates' names
                              epsilon  # equivalence limits
) {
  treatment.vec <- data.tested[, treatment]  # vector of the treatment variable
  treatment.group <- data.tested[treatment.vec == 1, ]  # treatment condition's data frame
  control.group <- data.tested[treatment.vec == 0, ]  # control condition's data frame
  results <- matrix(NA, length(covariates), 2)
  for (i in 1:length(covariates)) {
    results[i, ] <- equivarence.test(first.group = treatment.group, second.group = control.group, 
                                     covariate = covariates[i], epsilon = epsilon)
  }
  # p-value adjustment by the Benjamini-Hochberg method
  results[, 2] <- p.adjust(results[, 2], method = "BH")
  rownames(results) <- covariates
  colnames(results) <- c("std.dif", "p.value")
  results
}

# function to conduct F-tests for the analysis of primed labels
label.F.test <- function(variable, data, choice) {
  if (choice == TRUE) {
    lm.formula <- as.formula(paste0("choice ~ condition + ", variable))
  } else {
    lm.formula <- as.formula(paste0("rating ~ condition + ", variable))
  }
  lm.model.1 <- lm(lm.formula, data = data)
  if (choice == TRUE) {
    lm.formula <- as.formula(paste0("choice ~ condition * ", variable))
  } else {
    lm.formula <- as.formula(paste0("rating ~ condition * ", variable))
  }
  lm.model.2 <- lm(lm.formula, data = data)
  anova.result <- anova(lm.model.1, lm.model.2)
  c(anova.result$"F"[2], anova.result$"Pr(>F)"[2])
}

# function to conduct F-tests for the analysis of knowledge
knowledge.F.test <- function(variable, data, choice) {
  if (choice == TRUE) {
    lm.formula <- as.formula(paste0("choice ~ high.know + ", variable))
  } else {
    lm.formula <- as.formula(paste0("rating ~ high.know + ", variable))
  }
  lm.model.1 <- lm(lm.formula, data = data)
  if (choice == TRUE) {
    lm.formula <- as.formula(paste0("choice ~ high.know * ", variable))
  } else {
    lm.formula <- as.formula(paste0("rating ~ high.know * ", variable))
  }
  lm.model.2 <- lm(lm.formula, data = data)
  anova.result <- anova(lm.model.1, lm.model.2)
  c(anova.result$"F"[2], anova.result$"Pr(>F)"[2])
}

# function to compute marginal means controlling some variable
MM.control <- function(data, attribute, outcome, control, id) {
  lm.data <- data.frame(id = data[, id], A = data[, attribute], 
                        Y = data[, outcome], X = data[, control])
  l <- nlevels(lm.data$A)
  result.matrix <- matrix(NA, l, 3)
  for (i in 1:l) {
    level.label <- levels(lm.data$A)[i]
    result.matrix[i, ] <- unlist(lm_robust(Y ~ X, data = lm.data, 
                                           subset = A == level.label, 
                                           clusters = id)[c(1, 6, 7)])[c(1, 3, 5)]
  }
  result.matrix
}

#### load data ####
# task-level dataset
task.data <- read.csv("task_data.csv")

# reorder the levels of attribute variables
position.label <- c("Agree", "Neither", "Disagree")
task.data$a.article9 <- factor(task.data$a.article9, levels = position.label)
task.data$b.defense <- factor(task.data$b.defense, levels = position.label)
task.data$c.revisionism <- factor(task.data$c.revisionism, levels = position.label)
task.data$d.women <- factor(task.data$d.women, levels = position.label)
task.data$e.gay <- factor(task.data$e.gay, levels = position.label)
task.data$f.immigrant <- factor(task.data$f.immigrant, levels = position.label)
task.data$g.growth <- factor(task.data$g.growth, levels = position.label)
task.data$h.tax <- factor(task.data$h.tax, levels = position.label)

# respondent-level dataset
respondent.data <- read.csv("respondent_data.csv")

#### analysis conditioned by primed labels ####
## balance check (Table A.9)
round(EQT.wrapper.label(data.tested = respondent.data, treatment = "condition", 
                        covariates = c("gender", "age", "low.edu", 
                                       "middle.edu", "high.edu", 
                                       "DID.ratio", "ideology", 
                                       "article9.attitude", "defense.attitude", 
                                       "revisionism.attitude", "women.attitude", 
                                       "gay.attitude", "immigrant.attitude", 
                                       "growth.attitude", "tax.attitude", 
                                       "knowledge"), 
                        epsilon = 0.36), 3)

## compute MMs
choice.result.left.condition <- choice.result.right.condition <- 
  rating.result.left.condition <- rating.result.right.condition <- matrix(NA, 24, 3)
for (i in 1:8) {
  choice.result.left.condition[(3 * (i - 1) + 1):(3 * i), ] <- 
    MM(subset(task.data, condition == 1), issue.var.name[i], 
       "choice", "respondent.id")
  choice.result.right.condition[(3 * (i - 1) + 1):(3 * i), ] <- 
    MM(subset(task.data, condition == 0), issue.var.name[i], 
       "choice", "respondent.id")
  rating.result.left.condition[(3 * (i - 1) + 1):(3 * i), ] <- 
    MM(subset(task.data, condition == 1), issue.var.name[i], 
       "rating", "respondent.id")
  rating.result.right.condition[(3 * (i - 1) + 1):(3 * i), ] <- 
    MM(subset(task.data, condition == 0), issue.var.name[i], 
       "rating", "respondent.id")
}

# Figure A.8
cairo_pdf("Figure_A8.pdf", width = 6.2, height = 4, pointsize = 8)
par(mar = c(3, 0, 1.5, 0), lwd = 0.5, xpd = TRUE)
plot(NULL, NULL, type = "n", bty = "n", xlim = c(0.23, 0.85), ylim = c(0, 33), 
     xlab = "", ylab = "", xaxt = "n", yaxt = "n")
segments(seq(0.4, 0.6, 0.05), 0, seq(0.4, 0.6, 0.05), 33, lwd = 0.4, col = "gray")
for (i in 1:8) {
  segments(0.395, 35 - 4 * i, 0.605, 35 - 4 * i, lty = 3, col = "gray")
  segments(0.395, 34 - 4 * i, 0.605, 34 - 4 * i, lty = 3, col = "gray")
  segments(0.395, 33 - 4 * i, 0.605, 33 - 4 * i, lty = 3, col = "gray")
  segments(choice.result.left.condition[(3 * (i - 1) + 1):(3 * i), 2], 
           c(35 - 4 * i, 34 - 4 * i, 33 - 4 * i) + 0.2, 
           choice.result.left.condition[(3 * (i - 1) + 1):(3 * i), 3], 
           c(35 - 4 * i, 34 - 4 * i, 33 - 4 * i) + 0.2, col = plot.col[i])
  points(choice.result.left.condition[(3 * (i - 1) + 1):(3 * i), 1], 
         c(35 - 4 * i, 34 - 4 * i, 33 - 4 * i) + 0.2, 
         pch = 19, col = plot.col[i])
  segments(choice.result.right.condition[(3 * (i - 1) + 1):(3 * i), 2], 
           c(35 - 4 * i, 34 - 4 * i, 33 - 4 * i) - 0.2, 
           choice.result.right.condition[(3 * (i - 1) + 1):(3 * i), 3], 
           c(35 - 4 * i, 34 - 4 * i, 33 - 4 * i) - 0.2, col = plot.col[i])
  points(choice.result.right.condition[(3 * (i - 1) + 1):(3 * i), 1], 
         c(35 - 4 * i, 34 - 4 * i, 33 - 4 * i) - 0.2, 
         pch = 21, col = plot.col[i], bg = "white")
}
axis(1, at = seq(0.4, 0.6, 0.05))
segments(seq(0.65, 0.85, 0.05), 0, seq(0.65, 0.85, 0.05), 33, lwd = 0.4, col = "gray")
for (i in 1:8) {
  segments(0.645, 35 - 4 * i, 0.855, 35 - 4 * i, lty = 3, col = "gray")
  segments(0.645, 34 - 4 * i, 0.855, 34 - 4 * i, lty = 3, col = "gray")
  segments(0.645, 33 - 4 * i, 0.855, 33 - 4 * i, lty = 3, col = "gray")
  segments(0.75 + (rating.result.left.condition[(3 * (i - 1) + 1):(3 * i), 2] - 5.5) / 4, 
           c(35 - 4 * i, 34 - 4 * i, 33 - 4 * i) + 0.2, 
           0.75 + (rating.result.left.condition[(3 * (i - 1) + 1):(3 * i), 3] - 5.5) / 4, 
           c(35 - 4 * i, 34 - 4 * i, 33 - 4 * i) + 0.2, col = plot.col[i])
  points(0.75 + (rating.result.left.condition[(3 * (i - 1) + 1):(3 * i), 1] - 5.5) / 4, 
         c(35 - 4 * i, 34 - 4 * i, 33 - 4 * i) + 0.2, 
         pch = 19, col = plot.col[i])
  segments(0.75 + (rating.result.right.condition[(3 * (i - 1) + 1):(3 * i), 2] - 5.5) / 4, 
           c(35 - 4 * i, 34 - 4 * i, 33 - 4 * i) - 0.2, 
           0.75 + (rating.result.right.condition[(3 * (i - 1) + 1):(3 * i), 3] - 5.5) / 4, 
           c(35 - 4 * i, 34 - 4 * i, 33 - 4 * i) - 0.2, col = plot.col[i])
  points(0.75 + (rating.result.right.condition[(3 * (i - 1) + 1):(3 * i), 1] - 5.5) / 4, 
         c(35 - 4 * i, 34 - 4 * i, 33 - 4 * i) - 0.2, 
         pch = 21, col = plot.col[i], bg = "white")
}
axis(1, at = seq(0.65, 0.85, 0.05), labels = seq(5.1, 5.9, 0.2))
for (i in 1:8) {
  text(0.23, 36 - 4 * i, labels = issue.label[i], pos = 4)
  text(0.33, 35 - 4 * i, labels = "Agree", cex = 0.9, pos = 4)
  text(0.33, 34 - 4 * i, labels = "Neither", cex = 0.9, pos = 4)
  text(0.33, 33 - 4 * i, labels = "Disagree", cex = 0.9, pos = 4)
}
legend(0.23, 0, c("Left condition", "Right condition"), 
       pch = c(19, 21), bg = c(NA, "white"), bty = "n")
mtext(c("Choice", "Rating"), at = c(0.5, 0.75), cex = 1.2, font = 2)
dev.off()

## conduct F-tests
label.F.test.p.choice <- label.F.test.p.rating <- matrix(NA, 8, 2)
for (i in 1:8) {
  label.F.test.p.choice[i, ] <- label.F.test(issue.var.name[i], 
                                             data = task.data, choice = TRUE)
  label.F.test.p.rating[i, ] <- label.F.test(issue.var.name[i], 
                                             data = task.data, choice = FALSE)
}
rownames(label.F.test.p.choice) <- 
  rownames(label.F.test.p.rating) <- issue.var.name

# Table A.10
round(label.F.test.p.choice, 3)
round(label.F.test.p.rating, 3)

#### analysis conditioned by knowledge ####
## compute MMs
choice.result.high.knowledge <- choice.result.low.knowledge <- 
  rating.result.high.knowledge <- rating.result.low.knowledge <- matrix(NA, 24, 3)
for (i in 1:8) {
  choice.result.high.knowledge[(3 * (i - 1) + 1):(3 * i), ] <- 
    MM(subset(task.data, high.know == 1), issue.var.name[i], 
       "choice", "respondent.id")
  choice.result.low.knowledge[(3 * (i - 1) + 1):(3 * i), ] <- 
    MM(subset(task.data, high.know == 0), issue.var.name[i], 
       "choice", "respondent.id")
  rating.result.high.knowledge[(3 * (i - 1) + 1):(3 * i), ] <- 
    MM(subset(task.data, high.know == 1), issue.var.name[i], 
       "rating", "respondent.id")
  rating.result.low.knowledge[(3 * (i - 1) + 1):(3 * i), ] <- 
    MM(subset(task.data, high.know == 0), issue.var.name[i], 
       "rating", "respondent.id")
}

# Figure A.9
cairo_pdf("Figure_A9.pdf", width = 6.2, height = 4, pointsize = 8)
par(mar = c(3, 0, 1.5, 0), lwd = 0.5, xpd = TRUE)
plot(NULL, NULL, type = "n", bty = "n", xlim = c(0.23, 0.85), ylim = c(0, 33), 
     xlab = "", ylab = "", xaxt = "n", yaxt = "n")
segments(seq(0.4, 0.6, 0.05), 0, seq(0.4, 0.6, 0.05), 33, lwd = 0.4, col = "gray")
for (i in 1:8) {
  segments(0.395, 35 - 4 * i, 0.605, 35 - 4 * i, lty = 3, col = "gray")
  segments(0.395, 34 - 4 * i, 0.605, 34 - 4 * i, lty = 3, col = "gray")
  segments(0.395, 33 - 4 * i, 0.605, 33 - 4 * i, lty = 3, col = "gray")
  segments(choice.result.high.knowledge[(3 * (i - 1) + 1):(3 * i), 2], 
           c(35 - 4 * i, 34 - 4 * i, 33 - 4 * i) + 0.2, 
           choice.result.high.knowledge[(3 * (i - 1) + 1):(3 * i), 3], 
           c(35 - 4 * i, 34 - 4 * i, 33 - 4 * i) + 0.2, col = plot.col[i])
  points(choice.result.high.knowledge[(3 * (i - 1) + 1):(3 * i), 1], 
         c(35 - 4 * i, 34 - 4 * i, 33 - 4 * i) + 0.2, 
         pch = 19, col = plot.col[i])
  segments(choice.result.low.knowledge[(3 * (i - 1) + 1):(3 * i), 2], 
           c(35 - 4 * i, 34 - 4 * i, 33 - 4 * i) - 0.2, 
           choice.result.low.knowledge[(3 * (i - 1) + 1):(3 * i), 3], 
           c(35 - 4 * i, 34 - 4 * i, 33 - 4 * i) - 0.2, col = plot.col[i])
  points(choice.result.low.knowledge[(3 * (i - 1) + 1):(3 * i), 1], 
         c(35 - 4 * i, 34 - 4 * i, 33 - 4 * i) - 0.2, 
         pch = 21, col = plot.col[i], bg = "white")
}
axis(1, at = seq(0.4, 0.6, 0.05))
segments(seq(0.65, 0.85, 0.05), 0, seq(0.65, 0.85, 0.05), 33, lwd = 0.4, col = "gray")
for (i in 1:8) {
  segments(0.645, 35 - 4 * i, 0.855, 35 - 4 * i, lty = 3, col = "gray")
  segments(0.645, 34 - 4 * i, 0.855, 34 - 4 * i, lty = 3, col = "gray")
  segments(0.645, 33 - 4 * i, 0.855, 33 - 4 * i, lty = 3, col = "gray")
  segments(0.75 + (rating.result.high.knowledge[(3 * (i - 1) + 1):(3 * i), 2] - 5.5) / 4, 
           c(35 - 4 * i, 34 - 4 * i, 33 - 4 * i) + 0.2, 
           0.75 + (rating.result.high.knowledge[(3 * (i - 1) + 1):(3 * i), 3] - 5.5) / 4, 
           c(35 - 4 * i, 34 - 4 * i, 33 - 4 * i) + 0.2, col = plot.col[i])
  points(0.75 + (rating.result.high.knowledge[(3 * (i - 1) + 1):(3 * i), 1] - 5.5) / 4, 
         c(35 - 4 * i, 34 - 4 * i, 33 - 4 * i) + 0.2, 
         pch = 19, col = plot.col[i])
  segments(0.75 + (rating.result.low.knowledge[(3 * (i - 1) + 1):(3 * i), 2] - 5.5) / 4, 
           c(35 - 4 * i, 34 - 4 * i, 33 - 4 * i) - 0.2, 
           0.75 + (rating.result.low.knowledge[(3 * (i - 1) + 1):(3 * i), 3] - 5.5) / 4, 
           c(35 - 4 * i, 34 - 4 * i, 33 - 4 * i) - 0.2, col = plot.col[i])
  points(0.75 + (rating.result.low.knowledge[(3 * (i - 1) + 1):(3 * i), 1] - 5.5) / 4, 
         c(35 - 4 * i, 34 - 4 * i, 33 - 4 * i) - 0.2, 
         pch = 21, col = plot.col[i], bg = "white")
}
axis(1, at = seq(0.65, 0.85, 0.05), labels = seq(5.1, 5.9, 0.2))
for (i in 1:8) {
  text(0.23, 36 - 4 * i, labels = issue.label[i], pos = 4)
  text(0.33, 35 - 4 * i, labels = "Agree", cex = 0.9, pos = 4)
  text(0.33, 34 - 4 * i, labels = "Neither", cex = 0.9, pos = 4)
  text(0.33, 33 - 4 * i, labels = "Disagree", cex = 0.9, pos = 4)
}
legend(0.23, 0, c("High knowledge", "Low knowledge"), 
       pch = c(19, 21), bg = c(NA, "white"), bty = "n")
mtext(c("Choice", "Rating"), at = c(0.5, 0.75), cex = 1.2, font = 2)
dev.off()

## conduct F-tests
knowledge.F.test.p.choice <- knowledge.F.test.p.rating <- matrix(NA, 8, 2)
for (i in 1:8) {
  knowledge.F.test.p.choice[i, ] <- knowledge.F.test(issue.var.name[i], 
                                                     data = task.data, choice = TRUE)
  knowledge.F.test.p.rating[i, ] <- knowledge.F.test(issue.var.name[i], 
                                                     data = task.data, choice = FALSE)
}
rownames(knowledge.F.test.p.choice) <- 
  rownames(knowledge.F.test.p.rating) <- issue.var.name

# Table A.11
round(knowledge.F.test.p.choice, 3)
round(knowledge.F.test.p.rating, 3)

#### analysis considering candidates' display position ####
## effect of candidates' display position on outcomes
choice.result.displayed.column <- lm_robust(choice ~ position.right, 
                                            data = task.data, 
                                            clusters = respondent.id)
summary(choice.result.displayed.column)
rating.result.displayed.column <- lm_robust(rating ~ position.right, 
                                            data = task.data, 
                                            clusters = respondent.id)
summary(rating.result.displayed.column)

## compute MMs
# center the dummy variable of candidates' display position
task.data$position.right.centered <- scale(as.numeric(task.data$position.right), 
                                           scale = FALSE)

choice.result.column <- rating.result.column <- matrix(NA, 24, 3)
for (i in 1:8) {
  choice.result.column[(3 * (i - 1) + 1):(3 * i), ] <- MM.control(task.data, 
                                                                  issue.var.name[i], 
                                                                  "choice", 
                                                                  "position.right.centered", 
                                                                  "respondent.id")
  rating.result.column[(3 * (i - 1) + 1):(3 * i), ] <- MM.control(task.data, 
                                                                  issue.var.name[i], 
                                                                  "rating", 
                                                                  "position.right.centered", 
                                                                  "respondent.id")
}

# original results
choice.result <- rating.result <- matrix(NA, 24, 3)
for (i in 1:8) {
  choice.result[(3 * (i - 1) + 1):(3 * i), ] <- MM(task.data, 
                                                   issue.var.name[i], 
                                                   "choice", 
                                                   "respondent.id")
  rating.result[(3 * (i - 1) + 1):(3 * i), ] <- MM(task.data, 
                                                   issue.var.name[i], 
                                                   "rating", 
                                                   "respondent.id")
}

# Figure A.10
cairo_pdf("Figure_A10.pdf", width = 6.2, height = 4, pointsize = 8)
par(mar = c(3, 0, 1.5, 0), lwd = 0.5, xpd = TRUE)
plot(NULL, NULL, type = "n", bty = "n", xlim = c(0.23, 0.85), ylim = c(0, 33), 
     xlab = "", ylab = "", xaxt = "n", yaxt = "n")
segments(seq(0.4, 0.6, 0.05), 0, seq(0.4, 0.6, 0.05), 33, lwd = 0.4, col = "gray")
for (i in 1:8) {
  segments(0.395, 35 - 4 * i, 0.605, 35 - 4 * i, lty = 3, col = "gray")
  segments(0.395, 34 - 4 * i, 0.605, 34 - 4 * i, lty = 3, col = "gray")
  segments(0.395, 33 - 4 * i, 0.605, 33 - 4 * i, lty = 3, col = "gray")
  segments(choice.result.column[(3 * (i - 1) + 1):(3 * i), 2], 
           c(35 - 4 * i, 34 - 4 * i, 33 - 4 * i) + 0.2, 
           choice.result.column[(3 * (i - 1) + 1):(3 * i), 3], 
           c(35 - 4 * i, 34 - 4 * i, 33 - 4 * i) + 0.2, col = plot.col[i])
  points(choice.result.column[(3 * (i - 1) + 1):(3 * i), 1], 
         c(35 - 4 * i, 34 - 4 * i, 33 - 4 * i) + 0.2, 
         pch = 19, col = plot.col[i])
  segments(choice.result[(3 * (i - 1) + 1):(3 * i), 2], 
           c(35 - 4 * i, 34 - 4 * i, 33 - 4 * i) - 0.2, 
           choice.result[(3 * (i - 1) + 1):(3 * i), 3], 
           c(35 - 4 * i, 34 - 4 * i, 33 - 4 * i) - 0.2, col = plot.col[i])
  points(choice.result[(3 * (i - 1) + 1):(3 * i), 1], 
         c(35 - 4 * i, 34 - 4 * i, 33 - 4 * i) - 0.2, 
         pch = 21, col = plot.col[i], bg = "white")
}
axis(1, at = seq(0.4, 0.6, 0.05))
segments(seq(0.65, 0.85, 0.05), 0, seq(0.65, 0.85, 0.05), 33, lwd = 0.4, col = "gray")
for (i in 1:8) {
  segments(0.645, 35 - 4 * i, 0.855, 35 - 4 * i, lty = 3, col = "gray")
  segments(0.645, 34 - 4 * i, 0.855, 34 - 4 * i, lty = 3, col = "gray")
  segments(0.645, 33 - 4 * i, 0.855, 33 - 4 * i, lty = 3, col = "gray")
  segments(0.75 + (rating.result.column[(3 * (i - 1) + 1):(3 * i), 2] - 5.5) / 4, 
           c(35 - 4 * i, 34 - 4 * i, 33 - 4 * i) + 0.2, 
           0.75 + (rating.result.column[(3 * (i - 1) + 1):(3 * i), 3] - 5.5) / 4, 
           c(35 - 4 * i, 34 - 4 * i, 33 - 4 * i) + 0.2, col = plot.col[i])
  points(0.75 + (rating.result.column[(3 * (i - 1) + 1):(3 * i), 1] - 5.5) / 4, 
         c(35 - 4 * i, 34 - 4 * i, 33 - 4 * i) + 0.2, 
         pch = 19, col = plot.col[i])
  segments(0.75 + (rating.result[(3 * (i - 1) + 1):(3 * i), 2] - 5.5) / 4, 
           c(35 - 4 * i, 34 - 4 * i, 33 - 4 * i) - 0.2, 
           0.75 + (rating.result[(3 * (i - 1) + 1):(3 * i), 3] - 5.5) / 4, 
           c(35 - 4 * i, 34 - 4 * i, 33 - 4 * i) - 0.2, col = plot.col[i])
  points(0.75 + (rating.result[(3 * (i - 1) + 1):(3 * i), 1] - 5.5) / 4, 
         c(35 - 4 * i, 34 - 4 * i, 33 - 4 * i) - 0.2, 
         pch = 21, col = plot.col[i], bg = "white")
}
axis(1, at = seq(0.65, 0.85, 0.05), labels = seq(5.1, 5.9, 0.2))
for (i in 1:8) {
  text(0.23, 36 - 4 * i, labels = issue.label[i], pos = 4)
  text(0.33, 35 - 4 * i, labels = "Agree", cex = 0.9, pos = 4)
  text(0.33, 34 - 4 * i, labels = "Neither", cex = 0.9, pos = 4)
  text(0.33, 33 - 4 * i, labels = "Disagree", cex = 0.9, pos = 4)
}
legend(0.23, 0, c("w/ controlling position", "w/o controlling position"), 
       pch = c(19, 21), bg = c(NA, "white"), bty = "n")
mtext(c("Choice", "Rating"), at = c(0.5, 0.75), cex = 1.2, font = 2)
dev.off()